Localization threshold of Instantaneous Normal Modes from level-spacing statistics 
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We study the statistics of level-spacing of Instantaneous Normal Modes in a supercooled liquid. 
A detailed analysis allows to determine the mobility edge separating extended and localized modes 
in the negative tail of the density of states. We find that at temperature below the mode coupling 
temperature only a very small fraction of negative eigenmodes are localized. 
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I. INTRODUCTION 

The Instantaneous Normal Modes (INM) of a liquid are 
the eigenvectors of the Hessian (second derivative) ma- 
trix of the potential energy, evaluated at an instantaneous 
configuration. The interest in the equilibrium-average 
properties of the INM originates in the proposal to 
use them to study liquid dynamical properties, especially 
diffusion, which is considered to be linked to unstable 
modes (a subset of the modes with negative eigenvalue) 
JJ- They have been naturally apphed to ad- 
dress the problem of the glass transition: the glass phase 
is viewed as that where free (i.e. non- activated) diffusion 
is absent, and the disappearance of diffusion should be 
linked to that of the unstable modes Qj- The identifi- 
cation of these unstable modes presents some problems 

, and the localization properties of the INM are of 
interest. This is the problem we address here. 

Localization is an interesting and difficult problem in 
its own right. Given an x random matrix defined by 
the probability distribution of its elements in some (typ- 
ically site) base, the problem is to determine whether an 
eigenmode projects to an extensive number of base vec- 
tors (extended state) or not (localized state) in the large 
A^ limit. Only in a few cases there is a theoretical solu- 
tion for this problem [III IT^ . From the point of view 
of random matrix theory, the INM are the eigenvectors of 
an Euclidean random matrix (ERM) ^3 1 s-^d the prob- 
lem of localization in ERMs has been recently addressed 
analytically in the dilute limit 0| . Clearly, the problem 
is also hard from the numerical point of view, involving 
a question about the thermodynamic limit. Quantities 
such as the participation ratio |l5j , which can distinguish 
between localized and extended states but require knowl- 
edge of the eigenvectors, arc problematic numerically be- 
cause computation of eigenvectors is very expensive for 
large systems. 

Here we explore an approach [l^ relying only on eigen- 
values, based on the fact that the statistics of level spac- 
ing is strongly correlated with the nature of the eigen- 
modes (see sec. II). We apply it for the first time to a 
model glass-forming liquid at a temperature below the 
mode coupling temperature Tc |l7j . Our work can be re- 
garded as an extension of the results of ref.|^ as far as we 



perform a detailed analysis of the level spacing. Our em- 
phasis is on the exploring the usefulness of level-spacing 
statistics as a means to obtain a localization threshold 
in off-lattice systems and the possible limitations of this 
technique. 



II. THEORETICAL BACKGROUND 

The spectral function of the eigenvalues of a random 
matrix is S{\) = Tli^i^ ~ K)/N. Its (disorder-) aver- 
age is the density of states g{X) (DOS). The cumulative 
spectral function for this particular realization of the dis- 
order. 
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can be decomposed into a smooth part plus a fluctuating 
term 77fluct(A), whose average is zero. The smooth part is 
then 



C(A) EE (r?(A)) 



dX'g{X') 
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The level spacing is not studied directly on the Ai, 
because it depends on the mean level density. To elimi- 
nate this dependence, one "unfolds" the spectrum, which 
means to map the original sequence |Ai| onto a new one 
{Ci(Ai)} according to Eq. jSJ (see e.g. [13 for a detailed 
explanation of this point). The cumulative spectral func- 
tion can be expressed in terms of these new variables: 



,7(C)=r/(A(C)) =C + '7fluct(C)- 



(3) 



The distribution of the variable ( is uniform in the inter- 
val [0, 1] regardless the g{X). 

The nearest-neighbor spacing distribution P{s) gives 
the probability that two neighboring unfolded eigenval- 
ues and Ci+i are separated by s. It is one of the most 
commonly used observables in random matrix theory. It 
is different from the two level correlation function and it 
involves all the fc-level correlation functions with k > 2 
[T^ . It displays a high degree of universality, exhibiting 
common properties in systems with very different spec- 
tra. Although no general proof has been given, its shape 
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is thought to depend only on the locahzation properties 
of the states [l^l- In the case of the Gaussian Orthog- 
onal Ensemble (GOE), where all states are extended in 
the thermodynamic limit, it is known [TeI that P{s), nor- 
malized such that (s) = 1, follows the so-called Wigner 
surmise (also known as Wigner-Dyson statistics), namely 



WD 



is) 



TVS 



exp (-7rs^/4). 



(4) 



The linear behavior for small s is an expression of the 
level repulsion. This form actually characterizes many 
different systems with extended eigenstates (see e.g. 
ref . 12(1 and references therein). In the case of INM, it has 
been shown |^ that it describes the level spacing better 
and better as the fraction of localized states decreases. 

On the other hand, a system whose states are all local- 
ized will have completely uncorrelated eigenvalues. This 
corresponds to a Poisson process, and the statistics of 
two adjacent levels is given by 



Pp{s) = exp(-s). 



(5) 



If one deals with a set of levels which includes both local- 
ized and extended states, one expects some distribution 
interpolating between those two. A natural ansatz is the 
simple linear combination 



Plc{s; 7t) = {1- Tr)Pp{s) + ^Pwd(s), 



(6) 



which holds under the hypothesis that contributions com- 
ing from localized and extended modes simply add lin- 
early. Another possibility comes from a statistical argu- 
ment due to Wigner (see for example 0|), that leads to 
the heuristic function 



P{s) ~ fi{s) exp 



ds'fi{s') 



(7) 



where /i(s) is called level repulsion function. Taking 
fi{s) = CqS* with q G [0,1], one obtains the Brody dis- 
tribution 



PB(s;g) = Cqs'^exp 



r9+i[i/(<z + i)] 



q+l y " q+l 

(8) 

which interpolates between the Poisson {q = 0) and 
Wigner-Dyson {q — 1) distributions. However, this 
is just another phenomenological interpolation scheme, 
since there is no theoretical argument supporting a level 
repulsion function increasing as a power law with an ex- 
ponent smaller than one. 



III. METHOD 

The practical difficulty in performing the unfolding lies 
in finding a good approximation to the smooth (aver- 
aged) part of the cumulative spectral function, C(A). We 




FIG. 1: Evaluation of the level-spacing statistics. Top: INM 
spectrum of unit density soft-sphere (pair potential l/r^2) 
system at T = 0.68, as obtained from the numerical di- 
agonalization of 100 thermalized configurations. Middle: 
The cumulative function of the same system and decompo- 
sition in smooth and fluctuating parts (inset). Bottom: The 
level-spacing distribution of this system, normalized to have 
(s) = 1. Poisson and Wigner-Dyson distributions are also 
shown for comparison. 



have first obtained a cumulative function averaged over 
many samples of the Hessian (computed from a corre- 
sponding number of equilibrium configurations) and then 
taken (^(A) as the function defined by a cubic spline in- 
terpolation of the resulting staircase. Once this function 
is defined; the spacings of each sample can be evaluated 
by extracting the A values according to the g{X) and then 
computing s = (^(Aj+i) — C{\i); the histogram of these 
values is an estimate of the Pis) ■ We have also tried digi- 
tal filtering (Savistky-Golay [23 ) of the staircase, but the 
results were not satisfactory. The procedure is illustrated 
in Fig. HI 
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FIG. 2: Level spacing distributions Pi(s/(s)) and P2{s/{s)) 
obtained from the positive part of the spectrum of Fig. for 
several values of Ath- Wigner-Dyson (left) and Poisson (right) 
distributions are also plotted. 
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FIG. 3: DOS of GOE matrices at several iV. The solid line 
is the semicircular law predicted at A'' — > oo. Inset: zoom on 
the left tail. 



not affected by differences m sizes The problem of 
finding the threshold is then reduced to finding the maxi- 
mum oi Djs[Pi, P2] as function of Xth- We stress that the 
ideas behind the method arc justified only in the large 
limit, and a study of finite-size effects is thus crucial in 
this context. 



To estimate the localization threshold Al, we pro- 
ceed as follows. We divide the full spectrum into two 
parts at an arbitrary threshold Xth and study (after 
proper unfolding) the restricted level-spacing distribu- 
tions Pi(s) EE P{s\X < Xth) and P2(s) = P(s|A > A*,,). 
The localization threshold (where it exists) should corre- 
spond the value of Ath that leads to Pi{s) = Pwd(s) (ex- 
tended eigenstates) and ^2(5) = Pp{s) (locahzed eigen- 
states). On the other hand if Xth 7^ Al, Pi{s) and ^2(5) 
will bear more similarity to each other, since one of them 
will include spacings from both localized and extended 
levels. A qualitative feeling of what happens as Xth moves 
through the spectrum can be gathered from Fig. 13 where 
it can be clearly seen how P2is) evolves from a nearly 
Wigner to a Poisson distribution. We remark that these 
probabilities distributions are universal since no fitting 
parameters are required once the plot is versus s/{s), 
where (s) = JP{s)sds. 

Since one cannot say, based on a finite sample, when 
one of the distributions becomes "exactly" Poisson or 
Wigner, one looks for the value of Xth that makes both 
distributions "as different as possible from each other." 
To do this we use (following ref . IT^ the Jensen-Shannon 
(JS) divergence as a measure of the distance between two 
distributions. It is defined as 

D]s[Pi,P2]= H[aiPi+a2P2]-aiH[Pi]-a2H[P2]. (9) 

H[P] ~ — X^i log-P(-^j) is the Shannon entropy of 
the distribution P, and ai,a2 = 1 — ai are positive 
weights of each distribution. In what follows, we shall 
choose the weights as proportional to the support of the 
section of the (unfolded) spectrum considered to evaluate 
the level-spacing. This ensures that the JS divergence is 



IV. RESULTS 

A. A case study: the GOE 

We have applied this procedure to the GOE, as a test 
and illustration of the method. We generated ensem- 
bles of X iV random matrices for N = 10, 20, 50, 500 
with i.i.d. elements (taken from Gaussian distribution 
with zero mean and variance 1/^/N) and computed the 
DOS by numerical diagonalization (Fig. O. The JS di- 
vergence has a maximum that tends to the band edge as 
N grows (Fig.^ top), indicating that there is no localiza- 
tion threshold in this system, as it is known theoretically. 

To gain further insight into the workings of the 
method, we have also tried fitting the level-spacing dis- 
tribution restricted to eigenvalues lower than Xth with 
the functions interpolating between Poisson and Wigner- 
Dyson, thus defining a kind of order parameter for local- 
ization (tt in the case of the linear combination, Eq. |S1 
q in the case of the Brody distribution. Eq. IS)). Both tt 
and q should be zero if Xth ^ Al and non-zero otherwise. 
As Fig.^shows, both order parameters start from being 
different from their minimum at a value which roughly 
corresponds to the maximum of the JS divergence. How- 
ever, the minimum value is not zero, most likely due to 
finite-size problems. Unfortunately, it is not possible to 
verify that in this case, because increasing N decreases 
the fraction of localized states such that their number 
remains finite even when N ^ 00 Il9l . 
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FIG. 4: Top: JS divergence for GOE matrices at different 
A^. Middle: The Brody parameter (see text) for the same 
values of A'^. Bottom; The order parameter tt from the linear 
approximation ||SJ|. 
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FIG. 5: INM spectrum of the binary mixture of soft spheres 
at F = 1.49 as obtained from 300 equilibrium configurations 
(iV = 2048). 
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FIG. 6: The Jensen-Shannon divergence for the negative tail 
of the spectrum of Fig.|21 The DOS g{\) is also plotted (here 
it is normalized such that /^^5(A) — 1). 



B. The INM spectrum 

We have studied the soft-sphere binary mixture of 
ref. ^ at unit density and T = 0.2029 (to be com- 
pared with the mode-coupling critical temperature Tc ~ 
0.2262). Equilibration of the supercooled liquid at this 
temperature has been possible thanks to the fast Monte 
Carlo algorithm of ref. |2^ From the physical point of 
view, we are interested in studying the nature of negative 
modes (which represent about 4.3% of the total modes 
for this system). At the temperature considered, the dy- 
namics is highly arrested, and diffusion events are rare 
(indeed, at the mean field level, such as mode-coupling 
theory, diffusion is completely suppressed below Tc). Ac- 
cordingly, one expects that all or most of the negative 
modes correspond to localized eigenvectors (i.e. local re- 



arrangement of a non extensive number of particles). We 
find that this is not the case. 

In Fig. [5l we show the INM spectrum for a system of 
2048 soft spheres. The spectrum is expected to have 
two localization thresholds (on the positive and nega- 
tive tails), so to apply the scheme above we need first to 
separate the positive and negative modes. We focus on 
the negative modes. The JS divergence for the negative 
part, evaluated as explained above, is shown in Fig.Elfor 
N = 400,800 and 2048 particles. As N increases, the 
maximum of these curves does not shift as in the GOE 
example but it becomes sharper, pointing to a localiza- 
tion threshold. A quadratic fit of the peak, for the largest 
size system leads toAL = — 16.8±1.4. We also verify that 
the two distributions Pi{s) and ^2(5) are indeed Poisson 
and Wigner (respectively) for this threshold value. 

We next try to fit with the linear interpolation: in 
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FIG. 7: The parameter tv is found by fitting the level-spacing 
distribution P{s\X < \th) to the form in @. 



FIG. 8: The Brody order parameter for the supercooled liquid 
at different sizes. 



Fig. 13 we plot the fitting parameter tt (cf. Eq. 0) for the 
P{s\X < Xth)- As in the GOE case, the parameter goes to 
a non-zero value. At the values of N available to us, there 
is no clear evidence that larger system sizes will make tt 
go to zero for Xth < Al- To check whether this behavior 
is an indication of some non- linear effect, we have per- 
formed the following test. Assuming that the threshold 
is actually at Al, for each of the values of Xth of Fig. [7| 
we have generated random spacings distributed with the 
linear combination of Poisson and Wigner-Dyson. The 
weight TT was taken as proportional to the number of ac- 
tual levels between Al and Xth, i-e. 



TT (X I g{X)dX 



We then tried the same fitting procedure we applied to 
the INM spectrum, to see whether it would yield the 
same tt used. We found that samples of at least w 90000 
levels were needed in order to recover the correct weight 
with the fit (this is more than ten times greater than the 
number of levels available from the INM spectrum of the 
simulated liquid). Hence we attribute the finite value of 
TT (Fig. [TJ to finite-size effects. 

In the fit with the Brody distribution, the finite-size 
effects seem to be less pronounced (sec Fig. IS)). The 
Brody parameter is consistent with a localized phase for 
A < —16: here one can see that in the large N limit the 
order parameter goes to zero as A < Al. So the results 
from the fits and from the JS divergence are consistent 
with the existence of a localization threshold. 

Though a more accurate determination of the thresh- 
old needs larger system sizes, this results shows that most 
(more than 96%) of the negative modes in this system are 
of an extended nature. 



V. CONCLUSIONS 

The study of the level-spacing distribution of the INM 
spectrum of a glass forming liquid in the supercooled 
regime we have presented shows that it is possible to 
locate the mobility edge for the negative tail of this spec- 
trum with reasonable precision (other techniques like the 
inverse participation ratio are usually less precise). The 
INM level-spacing distribution is reasonably described in 
terms of Wigner and Poisson distributions and this in- 
formation can be used to determine the mobility edge. 

We have applied the technique to the soft-spheres bi- 
nary liquid below Tc- Our result can be summarized by 
stating that at this temperature only 3.4% of the nega- 
tive modes are localized. This adds to the evidence (see 
e.g. the critique by Gezelter et al. |1|) that not all ex- 
tended imaginary modes can be regarded as leading to 
free diffusion. It has been argued U 0, 0, Q that not all 
extended negative modes should be considered unstable 
in the sense of this approach [one should exclude false 
saddles (also called shoulder modes) and saddles that do 
not connect different minima, an analysis we have not 
done here] . Our result implies not only that many nega- 
tive modes cannot be regarded as diffusive (not surprising 
in view of earfier results, e.g. refs.0,|^, but also that the 
vast majority of these non diffusive negative modes are 
extended, even below Tc- 

It would be interesting to extend these results to study 
the temperature dependence of the localization proper- 
ties of the INM across the Mode Coupling temperature. 
We expect to do this in the near future. 
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